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Fine structural details of glycans attached to the conserved N-glycosylation site significantly not only affect function of 
individual immunoglobulin G (IgG) molecules but also mediate inflammation at the systemic level. By analyzing IgG 
glycosylation in 5,1 17 individuals from four European populations, we have revealed very complex patterns of changes 
in IgG glycosylation with age. Several IgG glycans (including FA2B, FA2G2, and FA2BG2) changed considerably with 
age and the combination of these three glycans can explain up to 58% of variance in chronological age, significantly more 
than other markers of biological age like telomere lengths. The remaining variance in these glycans strongly correlated 
with physiological parameters associated with biological age. Thus, IgG glycosylation appeal's to be closely linked with 
both chronological and biological ages. Considering the important role of IgG glycans in inflammation, and because the 
observed changes with age promote inflammation, changes in IgG glycosylation also seem to represent a factor contribut- 
ing to aging. 

Significance Statement Glycosylation is the key posttranslational mechanism that regulates function of immunoglob- 
ulins, with multiple systemic repercussions to the immune system. Our study of IgG glycosylation in 5,117 individuals 
from four European populations has revealed very extensive and complex changes in IgG glycosylation with age. The 
combined index composed of only three glycans explained up to 58% of variance in age, considerably more than other 
biomarkers of age like telomere lengths. The remaining variance in these glycans strongly correlated with physiological 
parameters associated with biological age; thus, IgG glycosylation appears to be closely linked with both chronological 
and biological ages. The ability to measure human biological aging using molecular profiling has practical applications 
for diverse fields such as disease prevention and treatment, or forensics. 
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AGING is a complex process of accumulation of molecu- 
lar, cellular, and organ damage, leading to loss of func- 
tion and increased vulnerability to disease and finally to death 
(1). It is well known that lifestyle choices such as smoking 
and physical activity can hasten or delay the aging process (2). 
Such observations have led to the search for molecular mark- 
ers of age that can be used to predict, monitor, and provide 
insight into age-associated physiological decline and disease. 



Protein structure is defined by the sequence of nucleo- 
tides in the corresponding genes, thus the polypeptide 
sequence of a protein cannot change with age. However, an 
important structural and functional element of the majority 
of proteins are the glycans that participate in virtually all 
physiological processes (3). Glycans are product of a com- 
plex pathway that involves hundreds of different proteins 
and are encoded in a complex dynamic network of hundreds 
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of genes (4). Epigenetic regulation of gene expression is 
expected to affect protein glycosylation and several publi- 
cations recently reported this effect (5-8). Changes in gly- 
cosylation with age have been shown over 20 years ago (9) 
and have also replicated in recent large population studies 
(10-13). 

Immunoglobulin G (IgG) is an excellent model glyco- 
protein because its glycosylation has been well defined 
(Figure 1), and many important functional effects of alterna- 
tive IgG glycosylation have been described (14). For example, 
glycosylation acts as a switch between pro- and anti-inflam- 
matory IgG functionality. Most of the IgG molecules are not 
sialylated and are proinflammatory. Terminal a2,6-sialylation 
of IgG glycans decreases the ability of IgG to bind to acti- 
vating FcyRs and promotes recognition by DC-SIGN, which 
increases expression of inhibitory FcyRIIB and is anti-inflam- 
matory (15). Another fascinating example is the role of core 
fucose in the modulation of antibody-dependent cellular cyto- 
toxicity: IgG-containing glycans that lack core fucose have 
100-fold increased affinity for FcyRIIIA and are therefore 
much more efficient in activating antibody-dependent cel- 
lular cytotoxicity than fucosylated glycoforms of the same 



molecule (16). On average, 95% of the IgG population is 
core fucosylated (12); thus, most of the immunoglobulins 
have a "safety switch," which prevents them from activating 
antibody-dependent cellular cytotoxicity. Malfunction of this 
system appears to be associated with autoimmune diseases 
as indicated by both pleiotropic effects of genes that asso- 
ciate with IgG glycosylation on different inflammatory and 
autoimmune diseases, and the observed alterations in IgG 
glycosylation in systemic lupus erythematous (17) and many 
inflammatory diseases (18). 

Interindividual variability of IgG glycosylation in a popu- 
lation is large (12) and it appears to be affected by both vari- 
ation in DNA sequence (19) and environmental factors (11). 
Most of the studies that investigated glycosylation changes 
with age were either of limited size or were performed on 
the total plasma glycome; thus, in addition to changes in 
glycosylation, the observed differences reflected changes 
in the concentration of individual plasma proteins. In this 
study, we focused on glycosylation of IgG and analyzed 
more than 5,000 individuals from four different European 
populations to provide definitive data about changes in IgG 
glycosylation through the lifetime. 
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Figure 1 . UPLC analysis of immunoglobulin G (IgG) glycosylation. Each IgG contains one conserved N-glycosylation site on Asn |t)7 of its heavy chain. Different 
glycans can be attached to this site and the process seems to be highly regulated. UPLC analysis can reveal composition of the glycome attached to a population 
of IgG molecules by separating total IgG N-glycome into 24 chromatographic glycan peaks (GP1-GP24), mostly corresponding to individual glycan structures. 
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Results 

Changes in IgG Glycosylation With Age in Four Large 
Human Populations 

IgG was isolated from plasma of 906 individuals from 
Croatian island Vis, 915 individuals from Croatian island 
Korcula, 2,035 individuals from Orkney Islands in Scotland, 
and 1,261 twins (560 monozygotic and 698 dizygotic) from 
the TwinsUK cohort. N-linked glycans were released, labeled 
with 2-aminobenzamide (2-AB), and analyzed by ultra per- 
formance liquid chromatography (UPLC) as described pre- 
viously (12). A representative chromatogram with structures 
of all identified glycans is shown in Figure 1. More detailed 
structural explanation is presented in Supplementary Table S 1 . 
After correction for multiple testing, statistically significant 
associations with age were observed for 21 out of 24 directly 
measured glycan structures (Table 1, Figure 2, Supplementary 
Figure SI). These associations replicated very well between 
the four populations, though some population-specific dif- 
ferences existed. The strongest association with age was 
observed in the level of galactosylation. Nongalactosylated 
glycans (A2 and FA2) steadily increased with age, whereas 
digalactosylated glycans (A2G2, FA2G2, A2BG2, and 
FA2BG2) decreased with age. However, with monogalacto- 
sylated structures (A2G1 FA2G1, A2BG1, and FA2BG1), the 
situation was more complex, with some glycans increasing 
and some decreasing depending on the position of galactose 
and the presence of bisecting GlcNAc (Supplementary Figure 
SI). The observed changes were quite different in men and 



women, with more dynamic changes in women (Figure 2). 
The most pronounced age associations in IgG glycosylation 
in women occurred between the ages of 45 and 55 when the 
majority of women enter menopause. All other elements of 
IgG glycosylation (fucosylation, sialylation, and bisecting 
GlcNAc) also strongly associated with age, sometimes in a 
quite complex manner (Supplementary Figure SI). 

To verify that longitudinal changes in IgG glycosylation 
within an individual follow trends observed at the popula- 
tion level, we have analyzed IgG glycosylation in 26 indi- 
viduals first sampled in year 2003 and then again in year 
2013. Two time points of relative levels of selected gly- 
cans in the same individuals are presented in Figure 3 and 
Supplementary Figure S2. In the majority of individuals, 
the change was very well aligned with changes observed at 
the population level, with women displaying more intensive 
changes than men. 

Prediction of Chronological Age From IgG 
Glycosylation 

Because most of the analyzed glycans strongly associ- 
ated with age, we attempted to build a predictive model of 
aging on the ORCADES cohort (2,035 individuals) using 
a multivariate analysis of covariance. The maximal model 
included 96 gender- specific glycomic parameters. The back- 
ward regression was performed and the best minimal model 
was determined according to the Akaike's information cri- 
terion (20). The minimal model was trained and validated 
by repeated random subsampling, using two thirds (1,357) 



Table 1 . Associations of Glycans With Age 



Glycan 




Korcula 




Orkney 




TwinsUK 




Vis 


R 


P 


R 


P 


R 


P 


R 


P 


GP1 


.128 


1.00E-04 


.372 


1.55E-57 


.264 


6.65E-18 


.024 


4.70E-01 


GP2 


.332 


7.31E-26 


.492 


3.13E-114 


.438 


1.59E-55 


.344 


7.64E-27 


GP4 


.555 


1.77E-85 


.645 


7.04E-266 


.622 


4.25E-152 


.574 


5.02E-87 


GP5 


.003 


1.73E-01 


.414 


1.07E-73 


.335 


5.69E-30 


.003 


9.29E-01 


GP6 


.553 


4.47E-81 


.713 


<lE-300 


.657 


6.84E-185 


.549 


5.43E-76 


GP7 


-.030 


2.35E-03 


.181 


3.68E-18 


.174 


1.09E-07 


.019 


4.37E-02 


GP8 


-.193 


6.49E-08 


-.171 


1.92E-12 


-.028 


6.86E-03 


-.231 


2.71E-11 


GP9 


.010 


3.91E-01 


.024 


3.89E-08 


.008 


2.3OE-03 


-.061 


9.06E-02 


GP10 


.141 


2.33E-07 


.335 


4.69E-46 


.183 


2.79E-10 


.043 


1.69E-01 


GP11 


.264 


1.05E-17 


.481 


4.31E-113 


.314 


3.77E-27 


.247 


5.36E-14 


GP12 


-.308 


2.44E-24 


-.227 


8.74E-36 


-.260 


2.21E-17 


-.241 


1.13E-14 


GP13 


-.151 


1.55E-07 


-.406 


2.60E-82 


-.277 


5.00E-21 


-.081 


1.59E-02 


GP14 


-.599 


9.20E-104 


-.716 


<lE-300 


-.656 


9.56E-182 


-.624 


6.75E-108 


GP15 


-.306 


1.55E-22 


-.354 


5.47E-64 


-.424 


6.67E-51 


-.432 


9.99E-42 


GP16 


.094 


5.82E-03 


-.005 


5.68E-06 


-.145 


1.37E-05 


.031 


1.55E-01 


GP17 


.03 1 


8.81E-02 


.005 


2.26E-09 


-.159 


3.75E-06 


.043 


4.55E-01 


GP18 


-.564 


2.86E-87 


-.665 


1.56E-295 


-.648 


6.40E-177 


-.607 


9.62E-100 


GP19 


.022 


5.07E-01 


.020 


4.05E-01 


.132 


1.49E-04 


.057 


9.11E-02 


GP20 


.005 


1.38E-01 


-.212 


6.35E-20 


-.121 


3.03E-03 


-.047 


2.27E-01 


GP21 


-.016 


1.14E-02 


-.067 


1.68E-02 


-.170 


3.49E-05 


.052 


1.25E-01 


GP22 


.033 


4.17E-02 


.033 


1.24E-02 


-.038 


2.05E-01 


-.002 


7.93E-02 


GP23 


-.313 


4.39E-21 


-.355 


3.83E-55 


-.293 


5.03E-21 


-.276 


4.99E-17 


GP24 


.032 


3.28E-01 


.114 


7.67E-06 


.058 


4.17E-02 


0029 


3.88E-01 



Note: GP1-24 - glycan peak 1-24; R = correlation coefficient; p = p-value. 
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Sex 



Figure 2. Relationship between age and selected glycan structures. Plots indicate associations between the individual contributions of four selected glycan 
structures to the total immunoglobulin G glycome and chronological age in four different human populations. Blue and red curves are fitted local regression models 
describing gender-specific relationship between age and glycan structure. The shaded region is a pointwise 95% confidence interval on the fitted values (there is 95% 
confidence that the true regression curve lies within the shaded region). 



of samples as the "training set" and one third (678) of sam- 
ples as the "validation set." In 1,000 replicas, regression 
coefficients were estimated on the training set and goodness 
of fit was obtained by applying the model estimated on the 
training set to the validation set. Final model estimates were 
calculated as median values of 1,000 iterations of random 
subsampling. The resulting GlycanAge index is composed 
of three glycan variables: 

Predicted male age = 139.9 + (85. 1 x GP6) - (5.2 x GP6 2 ) - 
(34.6 xGP14) + (l 1.8 xGP15) 

Predicted female age = 1 10. 1 + (164.5 x GP6) - 

(46.7 x GP6 2 ) - (22.4 x GP14) - 
(1.9xGP15) 



All three glycan structures selected by the optimal model 
were shown to be strongly associated with age, but each 
of the three glycan structures revealed different patterns 
between genders. The model explained 58% of variation in 
chronological age (with 95% of replicas giving R 2 between 
54% and 61%), with correlation between age and predicted 
age of .76 (.73-78) (Table 3). When applied to the full 
Orkney data set, the standard deviation (SD) of residuals 
was 9.7 years (Figure 4A). The predictive power of our 
GlycanAge index was better for women (64% of variance 
explained) than for men (49% of variance explained). 

The same model was tested on the remaining three popu- 
lations (Vis = 906, Korcula = 915, and TwinsUK = 1,261). 
Samples from these populations were processed in the 
same way as samples from Orkney and used to predict age 
with a model that was trained on the Orkney cohort. The 
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Figure 3. Individual changes in immunoglobulin G {IgG) glycans with age. Two time points of relative levels of glycans (percentages of individual glycan struc- 
tures within the total IgG glycome), which associated strongly with age within an individual, are presented on the background showing levels of the same glycans for 
different individuals from the same population. 
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25 50 75 25 50 75 

GlycanAge (years) GlycanBioAge (years) 

(A) (B) 

Figure 4. Prediction of chronological age from measured immunoglobulin G glycans. Using the simplified GlycanAge model, predicted age was calculated for 
2,035 individuals from the Orkney cohort and plotted against real chronological age of the same individual (A). The prediction of age was further improved by inclu- 
sion of forced expiratory volume and systolic blood pressure (B). 



new predictions based on the three independent popula- 
tions were somewhat lower, with correlation between age 
and predicted age of .64, .64, and .69 for Vis, Korcula, and 
TwinsUK, respectively. To see whether underperformance 
of prediction in these three populations was due to over- 
specificity of the trained model, starting from the same 
minimal model as in ORCADES and repeating the same 
procedure, we built additional models for each of the 
cohorts. Three models, trained separately on Vis, Korcula, 
and TwinsUK cohort, showed only marginally better pre- 
diction performance on their training cohorts (R = .65, .65, 
and .70) compared with the model trained on the Orkney 
cohort. This indicates two things. First, different popula- 
tions vary slightly in their associations between age and 
glycans. Glycans are more strongly associated with age in 
the Orkney population than in other three populations, and 
so accuracy of the age prediction in Orkney population is 
the highest. Secondly, although each population shows a 
different strength of association between age and glycans, 
glycan profiles change in similar way through lifetime 
across all four populations. As a consequence, a model 
trained on one population is capable of explaining almost 
all of existing association between age and glycans in 
another population. In our case, a model trained on Orkney 
cohort explained 42%, 41%, and 48% of variation of age in 
Vis, Korcula, and TwinsUK cohorts compared with 43%, 
43%, and 50% of variation explained by models trained 
specifically on these three populations. 

The same model was applied to a subpopulation of indi- 
viduals from Vis who were sampled again in 2013, 10 years 
after the original sampling in 2003. Although chronological 



difference between two samplings was 10 years for all 
individuals, the median value of age difference predicted 
from the GlycanAge index was 9.6 years for women and 
0.6 years for men (Supplementary Table S2). 

IgG Glycosylation and Biological Age 

To identify factors that may be responsible for the remain- 
ing variability in the GlycanAge index, we performed an 
association analysis with all available biochemical and 
physiological traits in our databases. Associations with sta- 
tistical significance after correction for multiple testing are 
shown in Table 2. Virtually all traits with strong association 
with GlycanAge in one or more studied populations (insu- 
lin, HbAlc, BMI, triglycerides, etc.) are known to be asso- 
ciated with unhealthy lifestyles. 

Both glycans and chronological age correlated significantly 
with a number of these parameters (Supplementary Table S3); 
therefore, we attempted to build a model that would combine 
biological information in glycans and other biological param- 
eters. The inclusion of forced expiratory volume in the first 
second (FEV1) and systolic blood pressure into the model sig- 
nificantly improved the prediction of chronological age. The 
extended model was trained and validated in the same way as 
the GlycanAge index and explained 71% (68%-74%) of vari- 
ation in chronological age of the Orkney cohort, with correla- 
tion between age and predicted age of .84 (.83-86). Just as 
for glycan age, the predictive power of this model was better 
for women (76% of variance explained) than for men (64% of 
variance explained). Using a minimal model constructed of 
two glycans and two biological parameters was 
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Table 2. Association of GlycanAge Index With Biochemical and 
Physiological Traits After Correcting for Chronological Age and Sex 



Orkney Vis and Korcula 





Beta 


p 


Beta 


p 


Insulin 


0.0755 


9.22E-08 


0.0402 


3.50E-01 


Fibrinogen 


0.0157 


1.98E-06 


0.0167 


8.83E-05 


HbAlc 


0.1106 


2.63E-06 




3 16E-03 


BMI 


0.0585 


1.67E-04 


0.0344 


1.04E-02 


Triglycerides 


0.0092 


1.75E-04 


0.0140 


1.20E-04 


Glucose 


0.0113 


2.09E-04 


0.0091 


4.77E-02 


Waist circumference 


0.1468 


2.08E-04 






Calcium 


0.0010 


2.35E-04 


0.0002 


7.04E-01 


D-dimer 


2.9670 


8.24E-04 






Cholesterol 


0.0036 


3.07E-01 


0.0201 


5.51E-08 


LDL 


0.0031 


3.26E-01 


0.0146 


6.08E-06 


Uric acid 


1.0773 


4.02E-02 


0.7620 


9.68E-04 


Note: HbAlc = £ 


;lycosylated hemoglobin; BMI 


= body 


mass index; 


LDL - low-density lipoprotein; p = 


p value; beta - regression coefficient. 


Predicted male a; 


ge = 49.6 


+ (9.8xGP6) 


-(2.2xGP14 2 )- 



(1 .2 x FEV1 2 ) + (0. 1 x SystolicBP) 

Predicted female age = 53.9 + (8.3 x GP6) - (2.4 x GP14 2 ) 
- (1 .9 x FEV1 2 ) + (0. 1 x SystolicBP) 

The model was tested on the Korcula cohort and the 
correlation between age and age predicted with the model 
trained on the Orkney cohort was .80. The model trained on 
the Korcula cohort explained 65% of variation in age in that 
cohort (70% in women and 57% in men), with a correlation 
of chronological and predicted age of .81 (.83 for women 
and .76 for men). An overview of all results, together with 
corresponding 95% range of replicas, can be seen in Table 3. 

Discussion 

Fine details of glycan structures attached to IgG Fc promote 
binding of IgG to different receptors and in this way modulate 
action of the immune system (14,21). Individual variability 
in IgG glycosylation is large (12) and this variability in gly- 
cosylation contributes to individual differences in function of 
the immune system (22,23). As it was clearly demonstrated 
for intravenous immunoglobulin therapy even slight changes 
in IgG glycosylation can direct pro- and anti-inflammatory 
actions of immunoglobulins (24). IgG Fc glycosylation mod- 
ulates inflammation (25,26) and through promotion or sup- 
pression of inflammation, it may significantly contribute to 
the process of biological aging (27-29). This hypothesis is 
further supported by the observed decrease in galactosylation 
in some premature ageing syndromes (13,30). 

Up to 50% of plasma glycome variability is estimated 
to be heritable (12,31). The remaining variability is appar- 
ently caused by environmental factors, but the observed 
long-term stability of the glycome (32) argues in favor of 



long-term epigenetic regulation of the glycan biosynthesis 
pathways (33). As we have shown in this study, a large part 
of this nongenetic variability can be explained by age and 
physiological variables related to age. In the total plasma 
glycome, as opposed to the IgG glycome, age was reported 
to have only minor effects on most of the glycans (31). 
Nearly all IgG glycans (21 out of 24 directly measured), 
however, were significantly affected by age (Table 1). The 
heavily influenced feature of IgG glycosylation was galac- 
tosylation, which decreased to less than 50% of its maximal 
value through lifetime. The decrease in IgG galactosyla- 
tion with age was initially reported more than 25 years ago 
(9) and replicated in a number of subsequent studies (34). 
Recently we reported that in childhood and adolescence, 
the galactosylation of IgG glycans actually increases with 
age (35). Indeed, from the population studies, it appears that 
total galactosylation is increasing until early adulthood (mid 
to late 20s) and decreases afterward. However, the situation 
is apparently more complex because the pattern of changes 
in men and women differ considerably (Figure 2). Also, 
galactosylated glycans with (GP15) and without bisecting 
GlcNAc (GP14) have different patterns of change. 

Galactosylation strongly decreases proinflammatory 
function of IgG (15); thus, both age-related and individual 
variation in the decrease of IgG galactosylation with age 
is affecting inflammation and can actively contribute to 
the deterioration of an aging organism through the inflam- 
mation-based process named inflammaging (27,36). The 
decrease in expression of the corresponding galactosyl- 
transferase ((34-GalTl) is clearly involved (37), but the rea- 
son and molecular mechanisms underlying such dramatic 
age changes in IgG galactosylation are not known. Our 
recent genome-wide association study of IgG glycosyla- 
tion identified polymorphisms in IL6ST, SMARCD3, HLA- 
DQA2/B2, and BACH2 that affect IgG galactosylation (17). 
None of these genes was previously reported to be involved 
in IgG glycosylation; thus, their functional role in this pro- 
cess is unknown. Recent study of genome-wide methyla- 
tion profiles identified a number of genes whose expression 
changed significantly with aging (38) and BACH2 was one 
of them. The combined observations that (i) BACH2 affects 
IgG galactosylation and (ii) that epigenetic regulation of 
BACH2 expression is changing with age suggest a new and 
intriguing hypothesis that changes in IgG galactosylation 
may actually represent a by-product of some other age- 
driven processes. In our data set, array methylation data for 
BACH2 and glycans were available for only 3 10 individuals, 
but even on this small data set and despite the "methylation 
noise" due to genetic and lifestyle factors, and DNA from 
other types of leukocytes, several significant associations 
between glycans and BACH2 promoter methylation were 
observed (Supplementary Table S4). With current knowl- 
edge, it is not possible to speculate whether associations of 
BACH2 promoter methylation changes with age and IgG 
glycosylation is causal or just a secondary by-product of 
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Table 3. Goodness-of-Fit and Spearman's Correlations of Chronological Age and Age Predicted by Various Models 




Training 


Test 


Population 


Female 




Male 




R 2 


Correlation 


R 1 


Correlation 


R 2 


Correlation 


GlycanAge Index 














Orkney 


Orkney 


57.8% [54.1-61.3%] 


.76 [.73 to .78] 


O^t.U 10 [_OU.U — Oo.l /C] 


on r 77 to 0 891 
.oU [_. / / LO U.oZJ 


AQ A<%-, \A\ Q Sfi 79M 


70 r f,A m 7S1 
. / U [.O^t LO . / J] 


Orkney 


Korcula 


41.3% 


.64 




. / J 




. JU 


Orkney 


Vis 


41.5% 


.64 


49.1% 


.70 


31.4% 


.56 


Orkney 


TwinsUK 


48.0% 


.69 


48.0% 


.69 


NA 


NA 


Korcula 


Korcula 


42.9% [34.2^19.7%] 


.65 [.58 to .70] 


51.1% [40.5-58.9%] 


.71 [.63 to 0.76] 


26.8% [16.1-36.3%] 


.51 [.40 to .60] 


Vis 


Vis 


43.0% [36.6-51.0%] 


.65 [.60 to .71] 


49.8% [41.3-57.8%] 


.70 [.64 to .76] 


34.6% [24.6-46.2%] 


.58 [.49 to .67] 


TwinsUK 


TwinsUK 


49.5% [42.9-55.3%] 


.70 [.65 to 0.74] 


49.5% [42.9-55.3%] 


.70 [.65 to .74] 


NA 


NA 


Combined Glycan Age Index 












Orkney 


Orkney 


71.3% [68.2-73.9%] 


.84 [.83 to 0.86] 


75.7% [72.7-78.4%] 


.87 [.85 to .89] 


63.7% [57.1-69.2%] 


.80 [.76 to .83] 


Orkney 


Korcula 


64.5% 


.80 


69.5% 


.83 


56.3% 


.75 


Korcula 


Korcula 


65.2% [59.7-70.0%] 


.81 [.77 to 0.84] 


69.6% [63.8-74.8%] 


.83 [.80 to .86] 


57.4% [46.5-67.0%] 


.76 [.68 to .82] 



Note: Numbers in brackets are 95% range of 1,000 replicas of random subsampling validation. 



BACH2 function. If changes in IgG glycosylation would be 
caused by changes in BACH2 promoter methylation status 
with aging, they would contribute to the deterioration of the 
organism with age by promoting inflammation, as proposed 
by Franceschi and colleagues (28,29). 

The accurate prediction of chronological and biologi- 
cal ages from biochemical parameters is a priority in the 
aging field (39-41). Here, we introduce a novel GlycanAge 
index that combines one nongalactosylated glycan (GP6) 
and two digalactosylated glycans (GP14 and GP15). This 
index predicts chronological age with error of 9.7 years and 
importantly explains 58% of variation in chronological age 
and sex. The explanation of nearly 60% of variance of age 
is impressive compared with other so-called age biomarkers 
like telomere length where most studies show 15%— 25% 
(42); thus, our new GlycanAge index appears to be more 
closely related to age than telomere lengths. After correct- 
ing for chronological age, the GlycanAge index correlated 
strongly with physiological parameters associated with 
measurements related to biological age (Table 2); thus, IgG 
glycosylation appears to be closely linked with both chron- 
ological and biological ages. We have further supported this 
by including forced expiratory volume in first second and 
systolic blood pressure in our prediction model, which sig- 
nificantly improved prediction and cumulatively explained 
71% of variation in chronological age (Figure 4B). 

The GlycanAge index, which includes only three gly- 
cans, seems to be a practical way to address both chron- 
ological and biological ages of an individual. The ability 
to measure human biological aging using molecular pro- 
filing has practical applications for diverse fields such as 
disease prevention and treatment, or forensics. The evalua- 
tion of aging immune system is here particularly interesting 
because glycosylation strongly affects function of immu- 
noglobulins (14) and the immune system in general (43). 
Whether aging causes changes in glycosylation of IgG or 
do glycosylation changes in IgG contribute to aging by pro- 
moting inflammation is very difficult to distinguish without 
longitudinal studies. Even though both of these aspects are 



probably partly true, it is clear that changes that occur in 
IgG glycosylation with aging are closely regulated and inte- 
grated with other physiological processes in the body. 

Materials and Methods 

Human Samples 

This study was based on plasma samples obtained from 
906 individuals (377 men and 529 women) from Croatian 
island Vis, 915 individuals (320 men and 595 women) from 
Croatian island Korcula (44), 2,035 individuals (797 men and 
1238 women) from the Orkney Islands in Scotland, and 1,261 
female twins (560 monozygotic and 698 dizygotic) from the 
TwinsUK cohort (45). The participants from Vis cohort were 
men and women between 18 and 88 years of age (median age 
56.0 and SD 15.2) and 18 and 93 years of age (median 56.5, 
SD 16.1), respectively. Male participants of Korcula cohort 
were aged between 20 and 90 (median age 57.8, SD 14.4) and 
female participants between 18 and 98 years of age (median 
55.4, SD 14.0). In Orkney cohort, men were aged between 
16 and 100 (median 53.7, SD 15.1) and women between 17 
and 97 years of age (median 52.9 and SD 15.3). TwinsUK 
cohort consisted of women aged between 27 and 83 years, 
with mean age of 58.6 and SD of 9.5. 

A subset of samples from the Croatian island Vis (n = 26) 
were sampled twice, first time in 2003 and second time in 
2013. All human participants included in this study signed 
an informed consent and the study was approved by the 
appropriate Ethics Committee of the University of Split 
Medical School, Institute for Anthropological Research, 
Research Ethics Committees in Orkney and Aberdeen, and 
Kings College London. This study conformed to ethical 
guidelines of the 1975 Declaration of Helsinki. 

Analysis of IgG glycans 

Isolation of IgG from human plasma. — The IgG was iso- 
lated using protein G monolithic plates (BIA Separations, 
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Ajdovscina, Slovenia) as described previously (12). Briefly, 
50-90 uL of plasma was diluted lOx with lx phosphate- 
buffered saline, pH 7.4, applied to the protein G plate and 
instantly washed with lx phosphate-buffered saline, pH 
7.4, to remove unbound proteins. IgGs were eluted with 
1 mL of 0. 1 M formic acid (Merck, Darmstadt, Germany) 
and neutralized with 1 M ammonium bicarbonate (Merck). 

Glycan release and labeling. — Glycan release and labe- 
ling of Korcula and Vis IgG samples was performed essen- 
tially as described by Royle and coworkers (46). Briefly, 
IgG was immobilized in a block of sodium dodecyl sulfate- 
polyacrylamide gel and N-glycans were released by diges- 
tion with PNGase F (ProZyme, Hay ward, CA). Each step 
was done in a 96-well microtiter plate to achieve the best 
throughput of sample preparation. After deglycosylation, 
N-glycans were labeled with 2-AB fluorescent dye. Excess 
of label was removed by solid-phase extraction using 
Whatman 3MM chromatography paper. Finally, glycans 
were eluted with water and stored at -20°C until usage. 

Orkney and TwinsUK IgG samples were first dena- 
tured with addition of 30 uE 1.33% sodium dodecyl sul- 
fate (w/v) (Invitrogen, Carlsbad, CA) and by incubation at 
65°C for lOmin. Subsequently, 10 uE of 4% Igepal-CA630 
(Sigma-Aldrich, St. Louis, MO) and 1.25 mU of PNGase 
F (ProZyme) in 10 uE 5x phosphate-buffered saline were 
added to the samples. The samples were incubated over- 
night at 37°C for N-glycan release. The released N-glycans 
were labeled with 2-AB. The labeling mixture was freshly 
prepared by dissolving 2-AB (Sigma-Aldrich) in dimethyl 
sulfoxide (Sigma-Aldrich) and glacial acetic acid (Merck) 
mixture (85:15, v/v) to a final concentration of 48mg/mL. 
A volume of 25 uE of labeling mixture was added to each 
N-glycan sample in the 96-well plate. Also, 25 uE of freshly 
prepared reducing agent solution (106.96mg/mL 2-picoline 
borane [Sigma-Aldrich] in dimethyl sulfoxide) was added 
and the plate was sealed using adhesive tape. Mixing was 
achieved by shaking for lOmin, followed by 2-hour incuba- 
tion at 65°C. Samples (in a volume of 100 uE) were brought 
to 80% acetonitrile (ACN) (v/v) by adding 400 uE of ACN 
(J.T Baker, Phillipsburg, NJ). Free label and reducing agent 
were removed from the samples using hydrophilic interac- 
tion chromatography-solid-phase extraction. An amount of 
200 uE of 0. 1 g/mL suspension of microcrystalline cellu- 
lose (Merck) in water was applied to each well of a 0.45 urn 
GHP filter plate (Pall Corporation, Ann Arbor, MI). Solvent 
was removed by application of vacuum using a vacuum 
manifold (Millipore Corporation, Billerica, MA). All wells 
were prewashed using 5 x 200 uE water, followed by equi- 
libration using 3x200 uE acetonitrile/water (80:20, v/v). 
The samples were loaded to the wells. The wells were sub- 
sequently washed seven times using 200 uE acetonitrile/ 
water (80:20, v/v). Glycans were eluted two times with 
100 uE of water and combined eluates were stored at -20°C 
until usage. 



Hydrophilic interaction chromatography-UPEC. — 
Fluorescently labeled N-glycans were separated by hydro- 
philic interaction chromatography on a Waters Acquity 
UPLC instrument (Waters, Milford, MA) consisting of a 
quaternary solvent manager, sample manager, and a FLR 
fluorescence detector set with excitation and emission 
wavelengths of 330 and 420 nm, respectively. The instru- 
ment was under the control of Empower 2 software, build 
2145 (Waters). Labeled N-glycans were separated on a 
Waters ethylene bridged hybrid (BEH) Glycan chroma- 
tography column, 100x2.1 mm i.d., 1.7 urn BEH parti- 
cles, with lOOmM ammonium formate, pH 4.4, as solvent 
A and acetonitrile as solvent B. Separation method used 
linear gradient of 75%-62% acetonitrile (v/v) at flow rate 
of 0.4mL/min in a 25-minute analytical run. Samples were 
maintained at 5°C before injection, and the separation tem- 
perature was 60°C. The system was calibrated using an 
external standard of hydrolyzed and 2-AB-labeled glucose 
oligomers from which the retention times for the individual 
glycans were converted to glucose units. Data processing 
was performed using an automatic processing method with 
a traditional integration algorithm after which each chroma- 
togram was manually corrected to maintain the same inter- 
vals of integration for all the samples. The chromatograms 
were all separated in the same manner into 24 peaks and the 
amount of glycans in each peak was expressed as percent- 
age of total integrated area. In addition to 24 directly meas- 
ured glycan structures, 53 derived traits were calculated 
as described previously (12). These derived traits average 
particular glycosylation features (galactosylation, fuco- 
sylation, and sialylation) across different individual glycan 
structures. Consequently, they are more closely related to 
individual enzymatic activities and underlying genetic pol- 
ymorphisms (47). 

Statistical Analysis 

To obtain normally distributed variables for 24 glycan 
structures, a log transformation was performed on glycan 
variables. Batch correction was performed for each glycan 
group separately using linear mixed model. In this model, 
dependent variable was log-transformed fraction of glycan 
group in total glycome. Age and gender were described 
as fixed effects, and technical sources of variation were 
described as random effects. To obtain corrected values for 
fractions of glycan groups, estimated batch effect (techni- 
cal source of variation) was subtracted from each log-trans- 
formed glycan measurement. 

The predictive model of age (GlycanAge model) was 
built using multivariate analysis of covariance approach 
implemented in "stats" package for R programming lan- 
guage (48). The maximal model included gender-specific 
linear and quadratic terms for each of 24 glycan struc- 
tures — 96 parameters in total. Feature selection was per- 
formed using backward elimination implemented in "step" 
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function of "stats" package for R, with Akaike information 
criterion (20) as optimization criteria. Reported goodness of 
fit, correlations, and regression coefficients were estimated 
as median values of 1,000 iterations of repeated random 
subsampling validation, with two thirds (1,357 samples) of 
Orkney cohort as "training set" and one third (678 samples) 
as "validation set." 

To identify factors that may be responsible for the dif- 
ferences between chronological age and predicted age, we 
performed association analysis with all available biochemi- 
cal and physiological traits in our databases. For each trait, 
we defined two models: 

predictedAge ~ age + sex 
predictedAge ~ age + sex + trait 

Likelihood ratio test was performed on those models to 
see if trait variable significantly improves model predic- 
tion. False discovery rate was controlled using Benjamini- 
Hochberg procedure. 

For combined biological and glycan age model, bio- 
chemical and physical parameters that were present in at 
least two cohorts were selected from available databases. 
The list contained many biomarkers that were previously 
shown to be related to human aging (49) and their corre- 
lations with age may be seen in Supplementary Table S3. 
Not all biomarkers were collected in all samples from 
Orkney, so resulting data set consisted of 1,728 people. The 
most significant biomarkers for the analysis were absent 
in Croatian Vis cohort and TwinsUK, so these populations 
were excluded from the analysis. 

The predictive model of combined biological and gly- 
can biomarkers was built in the same way as GlycanAge 
model. A full model consisted of gender-specific polyno- 
mial of second degree for each of the 24 biological param- 
eters included in the analysis, and three glycans included in 
GlycanAge, resulting in 108 parameters. Penalty for large 
number of parameters was selected based on 1,000 itera- 
tions of 10-fold cross-validation on Orkney's training set, 
with robustness of goodness of fit (adjusted for the num- 
ber of predictors) as selection criteria. The performance 
of the model was additionally tested on Korcula cohort. 
Training and validating on Korcula cohort were performed 
in the same manner. 

Supplementary Material 

Supplementary material can be found at: http://biomedgerontology. 
oxfordjournals.org/ 
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